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ABSTRACT 

Observed high-redshift QSOs, at z ~ 6, may reside in massive dark matter (DM) halos of more than 
10 12 Mq and are thus expected to be surrounded by overdense regions. In a series of 10 constrained 
simulations, we have tested the environment of such QSOs. The usage of Constrained Realizations has 
enabled us to address the issue of cosmic variance and to study the statistical properties of the QSO 
host halos. Comparing the computed overdensities with respect to the unconstrained simulations of 
regions empty of QSOs, assuming there is no bias between the DM and baryon distributions, and 
invoking an observationally-constraincd duty-cycle for Lyman Break Galaxies, we have obtained the 
galaxy count number for the QSO environm ent. We find that a clear discrepancy exists between the 
computed and observed galaxy counts in the lKim et al.l (|2009D samples. Our simul ations predi c t that 
on average eight z ~ 6 galaxies per QSO field should have been observed, while iKim et ail (2009) 
detect on average four galaxies per QSO field compared to an average of three galaxies in a control 
sample (GOODS fields). While we cannot rule out a small number statistics for the observed fields to 
high confidence, the discrepancy suggests that galaxy formation in the QSO neighborhood proceeds 
differently than in the field. We also find that QSO halos are the most massive of the simulated 
volume at z ~ 6 but this is no longer true at z ~ 3. This implies that QSO halos, even in the case 
they are the most massive ones at high redshifts, do not evolve into most massive galaxy clusters at 
z = 0. 

Subject headings: cosmology: dark matter — galaxies: evolution — galaxies: formation — galaxies: 
halos — galaxies: interactions — galaxies: kinematics and dynamics 



1. INTRODUCTION 

QSOs are among the most luminous objects in the uni- 
verse and can be detected at high redshifts, when the uni- 
verse was still very young — so fa r up to z ~ 6.42 (e.g., 
iFan et al.ll2003HWiliott et al.ll2010t ). The comoving space 
density at the bright end of the luminosity f unction (as 
detec ted by the Sloan Digital Survey, SDSS lYork et all 
20001) is - (2.2 ± 0.73) /i 3 Gpc~ 3 at z ~ 6 (jFan et al.l 
2004h . Observations of high-z QSOs raise a number of 
fundamental questions: where and how they formed, and 
what is their relationship to the formation of first stars 
and galaxies in the universe. High-z QSOs are important 
cosmological probes for studying the star formation his- 
tory, metal enrichment, early galaxy formation, growth 
of supermassive black holes (SBHs), properties of the 
interstellar and intergalactic matter, and the epoch of 
re-ionization in the universe. 

Rare, high-z QSOs have been claimed to form in highly 
overdense regions of the initial matter distribution. Their 
central SBHs are expected to have grown fast via mergers 
and/or accretion. Observed QSOs at z ~ 6 appear to 
host SBHs of ~ 10 9 M Q possibly residing within the most 
massive halos of a few xl0 12 M Q . Numerical simulations 
estimate a similar (to high-z QSOs) comoving density ~ 
1 Gpc~ 3 of these m assive halos fe.o.. lSpringel et aLll2005l : 
ISijacki et~a l. 2009). This similarity however, depends on 



the so-called QSO duty cycle, or a fraction of the time 
the QSO is actually active. If this duty cycle is less than 
unity, the QSOs are less rare and correspondingly reside 
in the less massive halos (e.g., Overzier et al. 2009). 
Additional arguments have been used in favor of high-z 
QSOs residing in lower mass halos (e.g., Willott et al. 
2 005). 

iKaiserl ()1984D and lEfstathiou fc Reesl (|1988t ) have 
shown that the presence of high peaks (rich clusters) 
in the primordial density field enhances their corre- 
lation function. Compared to clusters, galaxies have 
a smaller correlation length £(r) = 1 for r g ks 4 — 
7 ft -1 Mpc ()Davis fe Peebles! Il983f >. Therefore, primor- 
dial galaxy-size high peaks could in principle increase 
the clustering of galactic mass-scale halos in their vicin- 
ity (jMunoz fc Loebll2008h . 

Population of galaxies responsible for re-ionization of 
the Universe is not yet found. The epoch of re-ionization 
which extended from z ~ 15 has apparen tly ended by 
z ~ 6 , as indicated by observations (e.g., iBecker et al.1 
[20011 : iBarkanal l200l ICen fc McDonald! l2002h a nd by 
numerical modeling ( e.g.. iGnedin fc Ostrikerl 119971 : 
lHaiman fc H^^ded 120031 : iWvithe fc Loebll2003h . Because 
high overdensities will accelerate the galaxy formation 
and evolution, it is natural to look for these galaxies in 
the QSO neighborhoods. 
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A significant excess of sources compared to the density 
seen in GOODS 1 has been obtaine d from analyzing th e 
QSO J1030+0524 fi eld at z ~ 6.28 (IStiavelli et alJ[2005h . 
Zhen g et"aT1 (|2006T ) also observed a significant overden- 
sity 2 around the SDSS QSO J0836+0054 at z ~ 6. These 
observations suggest that galaxy clustering might win 
over a prospective negative feedback of QSO on its envi- 
ronment, and an excess of galaxies is associated with the 
QSO. However, a sample of Z775-dropout candidate galax- 
ies identified in five fields of the Hubble Space Telescope 
(HST) Advanced Camer a for Surveys (ACS) centered 
on SP SS QSOs at z ~ 6 (|Kim et al J 12003 : iMaselli et all 
2009) hint to a more complex behavior. Two fields has 
been claimed overdense, two underdense a nd an add i- 
tional one at the average density of GOODS (|Kim et al.D . 
IWillott et all (|2005t l detected' no ovcrdcnsities around 
SDSS QSOs, although their survey has been less sen- 
sitive than Kim's et al., including J1030+0524. 

In this paper, we model the formation of massive pure 
DM halos at z ~ 6 by means of Constrained and Uncon- 
strained numerical simulations (see section 2) and an- 
alyze the possible causes for the apparent discrepancy 
between simulations and observation s of hi gh-z QSO en- 
vironments described by iKim et ail (|2009T) . We choose 
~ 10 12 M & halos collapsing by z ~ 6 as both realistic 
and representative case. For this we resort to a sample 
of 10 different QSOs environments simulated at high- 
resolution. Our approach is different and yet comple- 
mentary to the study by Overzier et al. (2009), who have 
constructed the semi-analytical galaxy catalogues based 
on the Millennium simulation (Springel et al. 2005), and 
hence could make use of the sub-grid baryonic physics. 
The Millennium simulation volume 500 3 h _3 Mpc 3 is 
smaller than the typical volume occupied by a bright 
z ~ 6 QSO, if the QSO duty cycle is about unity (see 
above). However, it may provide a lower limit for what 
is expected for even more massive halos. 

Within the present cosmological formation scenario, 
ACDM, structure evolves in a 'hierarchical' way, with the 
growth of small density fluctuations amplified by gravity 
from an otherwise smooth density field. Within such sce- 
nario small objects collapse first and subsequently merge 
to form progressively larger and more massive structures. 
Therefore, the formation of high— 2 galaxies and QSOs 
depends on the abundance of DM halos within a given 
volume. 

Sim ilarity between the cosmic star form ation history 
(e.g., iMadau et aHl996t iBunker et al.l2004l) and the evo- 
lution of QSO abundances (e.g.. IShaver et al.lll996l ) sug- 
gests a possible link between galaxy formation and SBH 
growth. Evidence supporting this relation comes from 
the several correlations measured locally, i.e., at very 
low z, between the SBH masses and global properties 
of the host's sphe roid components, such as their masses 
and lu minosities (|Marconi fc HuntJl200l . light concen- 
tratio n dGraham et al.ll200lD andstellar velocity dispe r- 
sions (jFerrarese fc Merritdl20"00t iGebhardt et al.l l2000). 

The growth of SBHs can be plausibly linked to galaxy 

1 Great Observatory Origins Deep Survey (Giavalisco ct al. 

EH) 

2 As pointed out by Overzier (2009), this overdensity is less 
significant, due to the under-estimating the contamination from 
lower 2 interlopers. This comment applies also to Stiavelli et al. 
(2005) and Kim et al. (2009) 



formation proc ess. Metal-rich ga s associated with QSOs 
at high z (e.Q.. lBarth et al.ll2003| ) provides evidence that 
they are located at the centers of massive galaxies. It is 
likely, therefore, that high-z QSOs highlight the location 
of some of the first perturbatio ns that became nonlin- 
ear (e.o.. lTrenti fe Stiav elli 2007). Such objects reside in 
over-dense regions which might evolve into massive clus- 
ters of galaxies containing a population of large elliptical 
galaxies. At low redshifts, the S BHs of these elliptical 
galaxies might be dormant (e.g., iMagorrian et alj|1998t 
iLi et all [2007). If the SBHs evolve via gas accretion, 
their buildup can be regulated by a radiative and me- 
chanical feedback onto the infalling material. However, 
if the main mode of the SBH growth comes from mergers, 
the feedback is irrelevant. 

Observations so far appear to confirm the existence 
of some feedback effect on the QSO environment, but 
the results are not decisive. The large emission of radi- 
ation associated with the QSOs activity, might be suf- 
ficient to ionize the surrounding IGM and could even 
photoevaporate the gas from the neighboring DM ha- 
los, p reventing star formatio n before the gas cools down 
(e.g.. IShapiro fc Ra ga 2001) and suppressing galaxy for- 
mation in its vicinity. This could lead to a defi- 
ciency of (proto)galaxies around the QSOs despite the 
DM halo excess. On the other hand, the QSO ac- 
tivity could also lead to positive feedback, enhanc- 
ing the star for mation process and ther e fore, galaxy 
formation (e.g., |Begelman fc Cioffl 1 19891 iReesI 119891 : 
iShlosman fc Phinnevlll98l) ~ 

Finally, evidence for a positive feedback from the QSOs 
has been observed at lower re dshifts, e.g., for Ly man- 
break galaxies at z ~ 3 — 3.5 (jSteidel et al.l 120031 ) . and 
narrowband Ly q selected samples of rad io galaxies at 
2 < z < 5.2 (jVenemans et al.l |2003[ ). IStevens et all 
(2010) used mid-infrared imaging of five QSO fields and 
found evidence that the high-z QSOs are associated with 
a substantially elevated level of star formation activity 
at 1.7 < z < 2.8. Radio-quiet QSOs have shown a lower 
excess of star formation over radio-loud AGN at these 
redshifts. 

This paper is organized as follows: our methodology 
and numerics are explained in section 2, results are pro- 
vided in section 3, and comparison between observations 
and numerical simulations is given in section 4. Section 5 
discusses our results, and conclusions are given in the 
last section. The revised method for Constrained Real- 
izations implemented here is described in the Appendix. 

2. METHODOLOGY 

High-z QSOs are simulated by following a large cos- 
mological volume to accommodate the very low space 
density of such population. Such simulations exhibit a 
large dynamic ra nge to ensue the hie rarchical build up 
of the hosts (e.g., ISpringel et al.ll2005T ) and may include 
the gas dynam ics during mergers, star formation, etc. 
(ILi et all 120071) . 

An alternative approach is the use of the 
Constr ained Rea li zation s method (h er eafter 
CRs. iBertschingerl [l987t iHoffman fc Ribald U99H: 
Ivan de Wevgaert fc Bertschingerlll996| ) in order to pre- 
scribe the formation and collapse of a suitably massive 
DM halo at any redshift in a given computational box, 
with the subsequent addition of the baryonic component. 
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In this approach one can concentrate on the particular 
region of interest without the need to resort to large 
and expensive computati ons. Following th i s path , we 
adopt the CR method of iHoffman fc Ribald (|1991ft and 
impose the necessary constraints to seed DM halos of a 
fewxl0 12 M Q at z ~ 6, assuming that they harbor the 
observed high-z QSOs. 

2.1. Initial Conditions 

Our main goal is to assess the effects of a QSO-host 
DM halo on its immediate region at high— z. For this 
purpose we have created a suit of ten CRs with a DM 
halo seed of 10 12 Ii~ 1 Mq collapsing by z ~ 6, accord- 
ing to the top-hat model. As a control sample, we have 
constructed in tandem a suit of ten Unconstrained Real- 
izations (UCRs, hereafter). The same random seeds and 
cosmology have been used for each pair of CR - UCR re- 
alizations. Although ten realizations are not statistically 
significant, it is a large enough sample to overcome the 
cosmic variance. 

The CRs were constructed as follows. We have de- 
signed a set of ten different experiments, to probe 
different merging histories (i.e., environments) of a 
1O 12 /i _1 M0 halo in an ACDM cosmology. The cosmo- 
logical param eters are those from the 5 years WMAP 
data release (|Dunklev et al.l 12009ft . i.e., Q m = 0.279, 
fi A = 0.721, and h = 0.701, where h is the Hubble con- 
stant in units of 100 km s _1 Mpc _1 . erg = 0.817, the 
rms linear mass fluctuation within a sphere of radius 
8h~ 1 Mpc extrapolated to z = 0, was used to normal- 
ize the initial linear power spectrum. The constraints 
were imposed onto a grid of 256 3 within a cubic box of 
20 /i _1 Mpc. In the appendix |A1 we provide details on 
how the constraints have been imposed. All constraints 
were imposed at the center of the computational box. 

Figure [1] shows the typical outcome of the CR pro- 
cedure for one of the models. Each panel shows a 2-D 
cut (perpendicular to each other) along the central slice 
of the density field convolved with a Gaussian filter of 
10 12 Ii~ 1 Mq. The continuous lines represent overdense 
regions, the dashed lines — the underdense ones, while 
the thick continuous line is the average density within 
the box. The main constraint can be clearly noticed as 
the central peak in each panel. It is noteworthy to real- 
ize the presence of other smaller peaks surrounding the 
constrained region. These peaks, which are due to the 
random component of the CR method, could plausibly 
lead to major mergers at later epochs with the central 
overdensity. Smooth accretion of the surrounding mat- 
ter, together with mergers, make the halos grow con- 
tinuously, by z — 3 they have a mass ~ 10 13 M© (see 
Sec l3.2j) . Such halos, together with their surrounding en- 
viroir ment could potentially be the protocluster regions 
(e.g.. IVenemans et al.ll2007r i (see also Figure[3]). 

The UCR set was constructed using the same cosmo- 
logical parameters, grid and box sizes and random seeds, 
without imposing any constraint within the fields, result- 
ing in a truly corresponding unconstrained field of the CR 
one. Figure [2] depicts the central slice of a CR model (left 
panel) and its UCR counterpart (right panel) from our 
set. The normalization and contours for both fields are 
the same in order to stress similarities and differences be- 
tween the two models. Clearly absent is the peak at the 
center. However, the overall filamentary structure and 



underdense regions remain basically the same in both 
models. 

2.2. Numerical simulations 

The 10 + 10 simulations were performed with vacuum 
boundary conditions and physical coordinates. Due to 
these choices, a sphere of radius 10 /i _1 Mpc was carved 
out from each of the initial fields and evolved from 
z = 199 until z = 3 by means of the FT M-4.5 code 
([Heller fc Shlosmanl H99l iHeller et al] 12007ft . The mass 
resolution within our simulations, under the cosmologi- 
cal model assumed, is of m = 2.95 x 10 8 M Q . Therefore, 
a galactic mass halo of lO 11 ft _1 M0 or above can be re- 
solved with at least 1,000 particles. 

Figure [3] depicts the typical outcome for two of our sim- 
ulations, one CR (left panel) and its corresponding UCR 
realization (right panel) at z ~ 6. Both models corre- 
spond to the initial conditions maps shown in Fig. [2j 
Colors in both plots are proportional to their local par- 
ticle density and they have been equally normalized for 
a better comparison. The box size in both cases is 20 
Mpc. The most striking feature in the CR frame is the 
clustering of several halos of similar sizes and masses (see 
sec l3.1j) around the central peak. Such clustering is ab- 
sent in the UCR. One could naively assume that since 
there is only one constraint imposed onto the CR model, 
there should be only one halo of such mass present in 
the field. However, this is not the case as depicted in 
Fig. [31 a situation present in all of our models. The pres- 
ence of other equally massive or slightly smaller struc- 
tures around the QSO-constrained halo are due to the 
fact that the imposed constraint excites the formation 
of structures around it. This provides a higher density 
environment with respect to a "normal" environment, as 
shown in the right panel of Fig. |3l Such a behavior is 
a natural consequence of the hierarchical nature of the 
underlying CDM model. 

2.3. Halos 

DM hal os were identified by means of the HOP al- 
gorithm ()Eisenstein fc Hut] ll 998). This method isolates 
structures according to a purely density particle criteria. 
Densities have been calculated locally using an Smooth 
Particle Hydrodynamics (SPH) kernel. A structure or 
halo was defined as that region with a 3-D overden- 
sity contour 200 times the mean density. The lower cut 
in halo mass resolution was at N — 100 particles or 
2.95 x 10 10 M W to ensure tha t the derived halo masses 
are robust (|Trenti et al.ll2010l ). Figure [4] shows (central 
region of 6 Mpc) the typical outcome for the CR set. 
The circles are proportional to their virial mass. Notice 
that apart from the central halo with a mass > 1O 12 M , 
there are two other halos of the same mass order. Fur- 
thermore, the presence of several 10 11 and in particular 
of 1O 1O M0 halos is overwhelming. The latter halos fol- 
low the distribution of the more massive halos, indicating 
that the mass growth of such halos will continue substan- 
tially. The general excess of halos at any mass scale with 
respect to a "normal" (unconstrained) field can be easily 
noticed in their respective mass functions. 
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Fig. 1. — Initial density field for one of the CR models. Each panel depicts a perpendicular cut of a density cut along the central slice. 
Continuous lines represent overdensity regions, dashed lines — underdense regions, while the thick line is the mean overdensity within the 
box. The main constraint of 10 12 h~ Mq can be noticed at the center of each cut surrounded by the typical large scale structure of the 
Universe. The gray-scale vertical bar provides the linear density scale. 
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Fig. 2. — 2-D cuts of the central slice of thickness 78 h~ 1 kpc, of 
the initial density fields for one of our CR models (left panel) and 
for its respective UCR counterpart (right panel). Despite the main 
feature missing in the UCR model, the overall LSS is the same in 
both models. The contour lines arc the same for both models. 
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Fig. 4. — Halo distribution at z ~ 6 for a CR model (same as 
in previous figures). Circle sizes are proportional to their masses. 
The box size is 6 Mpc. 
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Fig. 3. — Typical 7V-bodyoutcome at z ~ 6 for a CR model 
(left panel) and its respective UCR (right panel). These models 
are the evolved fields of Fig. [2] Colors are proportional to their 
local densities. The box size is 20 Mpc for both models. 

3. RESULTS 

3.1. Halo mass function 

We have constructed halo Mass Functions (MF) from 
the halo catalogs for each model set, CR and UCR re- 
spectively. Figure [5] shows the total MFs averaged for the 



two sets of 10 realizations, and the compariso n Extended 
Press-Schechter curve fEPS lBond et~all[l99l in a spher- 
ical subvolume of radius 6 Mpc h^ 1 around the QSO halo 
location for each corresponding QSO - UCR simulation 
pair in our 10 + lOset. Predictions for the mass function 
of the QSO halos runs have been obtained by taking into 
account the average linear overdensity inside this sphere 
introduced by imposing the constraint. The presence of 
such overdensity changes the effective cosmology {e.g., 
Goldberg & Vogclcy 2004). The shaded area represents 
the 68% confidence interval, derived from the ±1<t varia- 
tion in the linear overdensity distribution from the set of 
10 CR runs. Noteworthy is the effect of the constraints 
imposed on the field at all mass scales. The MFs show 
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Fig. 5. — Average Mass Functions for the two sample sets within 
a sphere of radius of 6 h —1 Mpc. The continuous lines are the ex- 
pected number counts derived from the Extended Press-Schechter 
(EPS, Bond et al. 1991) formalism, for QSOs (red line) and UCRs 
(green line). The shaded areas represent the 68% confidence inter- 
val for a single realization of a CR run as derived from the linear 
overdensity variance among the sample of CR runs. Effectively, be- 
cause we are showing a mass function averaged over 10 realizations, 
the shaded area is more significant. 



that there is not only a difference at the high mass end 
range (imposed constraints regime), but the difference 
persists over the entire mass range. This result is in 
agreement with the hierarchical nature of the underlying 
CDM scenario. One would expect a better correspon- 
dence between the two MFs in the lower mass end, since 
the formation of smaller mass halos (M < 10 10 M Q ) is 
more common at these high redshifts. Indeed, such an 
agreement can be better noticed in Figure |6l where we 
split the MFs into two regions. The left panel represents 
the MFs of the inner box of 6 Mpc side centered around 
the main constraint (in the UCR cases we used the CR 
counterpart centers). Within this region the differences 
between both curves should be more striking. Outside 
this region, the corresponding MFs should have a better 
agreement as can be noticed on the right panel of Fig. [6l 
The two MFs differ in less than 5% overall (lower panel), 
which can be attributed to the fact that the region influ- 
ence by the constraint is somewhat larger than a radius 
of 3 Mpc. 

The corresponding inner MFs show that the environ- 
ment around a prospective QSO-host halo should be 
very rich in (sub)structure, at mass levels smaller than 
10 12 M Q . A massive halo appears not only to induce the 
formation of much smaller-scale, M < 10 10 Af Q halos, 
which one would expect to find here anyway, but also 
on all intermediate mass-scales which are rare at such 
redshifts. The relatively high amount of lO 11 ^^© halos 



(~ 10 on average) exceeds the expectations of those from 
unconstrained simulations, as shown by the total UCR 
MF (Figs. 02 EO right panel) . 

Figures [5] and [6] show that the potential QSO fields 
reside in anomalously dense environments. 

3.2. MAH: Fate of the QSO-host halos 

The next logical step is to ask what is the fate of such 
high-density regions? Do they evolve into the most mas- 
sive clusters nowadays or to a more "normal" cluster- 
group size halos? Majority of the performed work in 
this field indeed assumes that high-redshift QSOs can 
exist only in the precursors of what is now the most 
massive clusters (m > 10 15 M W ) in the universe (e.g., 
ISpringel et al.ll2005t iLi et al.ll2007j ). Such an assumption 
seems natural since the DM halos can only grow with 
time. 

However, it has been shown that this rule might 
not necessarily be true, both by means of the 
EPS formalism an d by the ana l ysis o f the Millen- 
nium si mulation (ISpringel et all |2005[) merge r tree 
history rtPe Lucia fc BlaizotJ l2007t iTrenti etaLl 120081 : 
lOverzier et al.ll2009D . The present DM mass of the QSO- 
host halos identified at z ~ 6 could be of the order of 
m ~ 10 14 Mq, and most massive halos identified at var- 
ious redshifts do not necessarily maintain this property 
at lower redshifts. 

This effect is confirmed in our simulations, although 
they have been stopped at z = 3. Throughout their 
evolution, the halos have experienced a series of ma- 
jor mergers and have grown to ~ 6 x 10 12 M Q . The 
most massive halo at z = 3, the QS09 has grown to 
~ 1O 13 M0, although it is far from being most massive at 
z ~ 6. We have constructed the Mass Accretion Histo- 
ries (MAHs) for each of the QSO-host halos in our CRs 
set. For this purpose, we have reproduced their respec- 
tive merger trees and followed the branch that represents 
the constrained halo, at z = 6, imposed by the initial 
conditions. 

Figure [7] represents the MAHs for all CR models. 
Most of them (apart from model QS09) follow the 
similar trajectories. The abrupt mass increases in the 
MAH curves are associated w ith major mergers [e.g., 
iRomano-Diaz et alJl20"o"l l2007j ). which seem to occur at 
all calculated redshift s. If we compare our M AHs with 
the cluster sample of iTasitsiomi et all ((2004) and their 
fitted MAH, we find that our halos will end up as a nor- 
mal cluster size halos. Furthermore, a simple estimate 
accounting for all the mass that could collapse into such 
DM halos from z = 3 till the present time, indicates that 
halos will grow to - 9 x 10 13 - 2 x 10 14 M Q . 

The natural dispersion of the MAH curves in Fig. [7] 
appears to be the result of the cosmic variance and comes 
from varying the seed of the initial conditions. A striking 
feature of these curves is their relative evolution. For 
example, the curve corresponding to QSO10 appears to 
be least massive of the sample at z ~ 6 while it grows 
fast and becomes the most massive at z ~ 3. 

4. COMPARISON WITH OBSERVATIONS: PREDICTIONS 
FOR GALAXY COUNTS IN QSO FIELDS 

Pure DM simulations can be used to predict the ex- 
pected Lyman B reak Galaxies number counts in the 
iKim et all (|2009[ ) observations by making additional as- 
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Fig. 6. — Average Mass Function for the two sample sets within a box of 6 Mpc (left panel), and outside this inner box (right panel). 
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factor a (lower axis). 



TABLE 1 



QSO Fields 



SI 



S2 



S3 



GOODS 
J1030+0524 
J1049+4637 
J1148+5251 
J1306+0356 
J1630+4012 
Mean 



8.08 ±2.84 
+5.92 ±4.70 
-0.08 ±4.02 
-5.08 ±4.02 
-7.08 ±4.02 
+2.92 ±4.37 
-3.40 ±1.89 



3.95 ± 1.99 
+4.05 ± 3.46 
-1.95 ±2.81 
-1.95 ± 2.81 
-3.95 ±2.81 
+4.05 ± 3.46 
+0.05 ± 1.38 



2.96 ± 1.72 
+7.04 + 3.60 
+1.04 + 2.64 
-2.96 + 3.60 
-1.96 + 2.43 
+2.04 + 2.82 
+1.04+1.37 



Note. — Number of £775 - dropouts and Pois son Error 
by S/N and Color Limit from IKim et"aTI l(2009h : (1) ob- 
served QSO fields; (2) extra-galaxies in the SI sample (i.e., 
the GOODS field subtracted from the actual galaxy count 
in Kim et al. [2009]); (3) extra-galaxies in the S2 sample; 
(4) extra-galaxies in the S3 sample. Except in the GOODS 
line, all numbers have GOODS fields subtracted. The SI 
— S3 samples are defined in the text. These numbers do 
not include the target QSO galaxies. The signs "+" mean 
overdensity and "— " mean underdensity with respect to the 
GOODS fields (first line). The 'Mean' value for S3 refers 
to mean galaxy counts in 5 observed fields (not to excess 
counts). 



sumptions. The QSO host galaxy is not counted in the 
following statistics. First, for simplicity, we assume that 
one galaxy is found per one DM halo. Second base d on 
the field of view geometry of the IKim et all (|2009| ) ob- 
servations, we compute the comoving volume probed in 
an ACS field of view, 11.3 arcmin 2 , i.e., 6x6 (h~ 1 Mpc) 2 
at z ~ 6, for i-dropouts. Assuming a pencil beam ge- 
ometry and redshift uncertainty of Az = 1 centered at 
z — 6, we obtain 6x6 x 317 (fe~ 1 Mpc) 3 fr o m the cosmic 
variance calculator 3 of iTrenti fc Stiavellil (|2008l) . The 
effective volume for the search is however smaller due 
to both incompleteness and because the faint i-dropout 
sources cannot be detected if they are found in the prox- 



imity (or behind) brighter foreground galaxies and stars. 
Based on artificial source recovery simulations carrie d 
out to search for z > 5 galaxies fe.o., IOesch et al.ll2007l) . 
we estimate that the effective volume is typically one 
half of the pencil be am volume. Hence, we assume that 
each ACS field in the lKim et al.l observations has probed 
~ 5,706 (h^Mpc) 3 . 

Two values for S/N color limits of 1775 — zgso = 1.3 
and 1.5 have been considered. Choice of S/N > 5 and 
*775 — z 850 > 1-3 have been labeled SI; S/N > 5 and 
*775 — z 850 > 1-5 have been labeled S2; and S/N > 8 and 
*775 — Z 85Q > 1-3 have been labeled S3. 



3 http: //casa. Colorado . edu/\lower3pt\hbox{\charl26}trenti/CosmicVariance .html 
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TABLE 2 



ID 


Duty-Cycle 




UCR 


CR (QSO) 


Extra-galaxies 


(Kim et al.) 


e DC 


(10 10 ft- 1 M ) 


Halos 


Halos 




S3 


1.0 


11.90 


0.15 


13.4 


13 


S2 


1.0 


11.10 


0.2 


14.5 


11 


SI 


1.0 


8.80 


0.6 


18.6 


18 


S3 


0.5 


10.60 


0.3 


15.2 


7 


S2 


0.5 


8.80 


0.4 


18.7 


9 


SI 


0.5 


6.38 


1.2 


28.1 


13 


S3 


0.2 


6.76 


0.7 


26.8 


5 


S2 


0.2 


5.84 


1.0 


31.2 


6 


SI 


0.2 


4.08 


3.0 


48.5 


9 


S3 


0.1 


4.66 


1.5 


40.9 


1 


S2 


0.1 


4.11 


2.0 


48.1 


5 


SI 


0.1 


2.83 


6.1 


69.7 


6 



Note. — Predicted galaxy counts: (1) the samples from lKim et alj (|2009T ); (2) the 
duty cycle; (3) minimum DM halo mass — M m i n ; (4) the UCR halos; (5) the QSO halos 
counts which represent the average number of halos in the respective sets more massive 
than M m in; (6) the extra number of galaxies expected. All values have been averaged 
over 10 CRs and 10 UCRs. 



Table Q] reproduces the o bservational counts in all 
QSO fields (|Kim et alJl2009t ) with one difference — the 
GOODS counts have been subtracted. So each number, 
with the exception of the GOODS field, represents the 
extra-galaxy counts. 

We use the galaxy counts in the GOODS field, as re- 
ported in the first data line of Table [1] to derive the 
number density of i-dropouts galaxies in a typical region 
of the universe at z ~ 6, depen ding on the different obser- 
vation selection procedures in iKim et all . These number 
densities have been compared to number densities of DM 
halos derived from our sample of unconstrained control 
simulations. 

The minimum DM halo mass, M m i n , required to match 
the observed number density of galaxies has been calcu- 
lated by assuming a given duty cycle 4 £dc for the Ly- 
man Break Galaxies. The resulting M m - m is reported 
in Table [2] for different values of cdc and for the dif- 
ferent selections SI - S3. The values of M m i n are con- 
sistent with those der ived from a conditional luminos- 
ity function appro ach (|Stark et al.ll2009t iLee et ai1l2009l ; 
iTrenti et al.ll2010l) . as well as from the c lustering prop- 
erties of i-dropouts (|Overzier et al.l [20061 ). 

Next, we have considered the CR runs and counted the 
average number of halos with M > M m i n in the prox- 
imity of a QSO host halo, using a pencil beam volume 
of 6 x 6 x 12 (/i^Mpc) 3 centered at the QSO position. 
The volume size is determined by the ACS field of view 
and by the limits imposed by the size of our simulation 
box. The depth of 12 h~ 1 Mpc guarantees that the vol- 
ume considered does not extend outside the simulation 
computational volume, which is progressively reduced as 
the overdensity imposed at the center of the box col- 
lapses. We have subtracted the number of UCR halos 
above M > M min (column 4 of Tabled from that of QSO 

4 The duty cycle is defined here as a probability, edCi ranging 
between and 1, for a given high-z galaxy to be a Lyman Break 
Galaxy 



halos in the pencil beam volume around the QSO (col- 
umn 5 of Table [2]) . After multiplying this excess number 
by cdc j we have obtained the excess number of galaxies 
predicted in a QSO field compared to the GOODS field 
(last column in Table [21 ■ For a typical duty cycle val- 
ues of 0.2 (e.g., [Wvithe fc Lqebll2006b IStark et all 120071: 
ILee et al.|[2009t ITrenti et al.ll2010t ). we expect ~ 5 - 9 
more galaxies in the QSO fields, depending on the ob- 
servational selection criteria, SI - S3. We note also that 
while the above value of cdc appears as a peak in the 
likelihood countq urs, it is a very shallow maximum (e.g., 
IStark et al J 120071) . A variation in edc value between 0.1 
and 0.5 does not substantially change the result of the 
galaxy counts within each of the samples, SI - S3. 

The galaxy counts in Table [5] do differ when various 
samples are compared. If we restrict ourselves to the 
conservative S3 sample, because it has the highest S/N 
value, the predicted overdensity corresponds to ^ 5 extra 
galaxies, averaged over 10 QSO fields. 

Comparison between the observational counts in Ta- 
ble [1] and predictions from numerical simulations in Ta- 
ble [D shows that a broad disagreement exists between 
both counts of extra galaxies. Observational counts in 
each sample are consistent with no galaxy overdensity or 
underdensity in the QSO fields, except perhaps for QSO 
J1030 field (overdensity) and J1148 (underdensity). If we 
put more weight on S3, because it has the highest S/N 
value, this conclusion strengthtens even more. The aver- 
age value of galaxy counts in QSO fields in excess of the 
GOODS fields for the S3 selection is 1.04 ± 1.37 which is 
much smaller than the predicted 5 galaxies. This clearly 
suggests that galaxy formation around a QSO halo does 
not follow the predictions based on the abundance of DM 
halos in this environment compared to the field. This dis- 
crepancy is further discussed in section 5. A similar con- 
clusion has been also reached by Overzier et al. (2009) 
based on the counts of simulated dropout galaxies. 

Finally, we note that model halos around the QSO are 
distributed asymmetrically due to the presence of fila- 
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ments converging at the QSO location (see Figure SJ, 
hence the asymmetry of sources observed in Fig. 6 of 
iKim et all is a natural consequence of the topology of 
DM halo distribution around the QSO halo. 

5. DISCUSSION 

Observations bv lKim et all (j2009h appear broadly con- 
sistent with there being no overdensity or underdensity 
in the galaxy counts of QSO fields with respect to the 
GOODS fields. This result is puzzling in view of the sub- 
stantial overdensity of DM halos in the vicinity of the 
QSO host halo in the numerical simulations presented 
here. Taken at face value, the observations by Kim et al. 
imply that, at a given mass, DM halos in the vicinity of 
a bright QSO host less luminous galaxies than DM halos 
with the same mass i n the contr ol GOODS fields. 

On the other hand, IKim et all find that the 1775 — zsso 
color distribution of dropout galaxies differs significantly 
from the averages for the GOODS galaxies similarly se- 
lected, at 99% confidence level. This suggest that de- 
spite having a comparable average number of galaxies, 
the QSO fields have observable differences compared to 
the GOODS fields. 

The simplest resolution of this discrepan cy may lie 
in the relatively low number statistics of the IKim et all 
(2009) observations. Only five QSO fields were observed 
and the number of galaxies detected per field is low (3- 
8 on average depending on the selection adopted). In 
fact, if the galaxy number counts for the QSO fields 
with S3 selection were drawn from a Poisson distribu- 
tion with < N >= 8 (3 galaxies from GOODS plus pre- 
dicted excess of 5 galaxies), then there is a probability 
p ~ 0.04 of measuring on average 4 galaxies per QSO 
field over the five fields observed. Cosmic variance is 
likely to increase the field to field variations, increasing p. 
But on the other hand, if the QSO halo i s more massive 
than 10 12 M(7) (as for example suggested by lSpringel et al.l 
120051 ). ~ 1O 13 M (e.g., Dijkstra et al. 2008), the excess 
number of galaxies in its surroundings will also be higher 
(see also Munoz & Loeb 2008). 

O f course, if the QSO halo is less massive instead, then 
the 1 im et al.l results are instead more likely. As we 
stated in the Introduction, a number of arguments do 
point to QSOs residing in the less massive halos (e.g., 
Walter at al. 2004; Willott et al. 2005; Overzier et al. 
2009). We note, that although the duty cycles of AGN, 
in general, and QSOs, in particular, are hotly debated 
issues at present, the longest duty cycle among AGN 
appears to be ~ 10 8 yrs — that of radio galaxies. This 
hints at the high-z QSO duty cycle being less than unity. 

It is possible that the projected (on the sky) overden- 
sities around the high-z QSOs are not associated with 
the actual 3-D overdensities, and are the result of con- 
tamination by the 'field' dropout galaxies (e.g., Overzier 
et al. [2009]). In their Fig. 15, Overzier et al. exhibit 
~ 30 projected overdensities from the Mock Survey of 
~ 70,000 arc min 2 area based on the Millennium simu- 
lation. Many of these overdensities do not correspond to 
the 3-D overdensities. Therefore, we estimate the prob- 
ability of such contamination in the Kim et al. (2009) 
fields. Kim et al. used the ACS fields which correspond 
to 1/30 of a 316 arc min 2 GOODS field. Hence, each 
Kim et al. field is ~ 316 x 30 = 9,480 times smaller 
than the Overzier et al. field. Assuming here the worst 



case scenario that all of the 30 overdense regions in their 
Fig. 15 are 'fake' projection cases, we obtain the prob- 
ability of such a fake overdensity in Kim et al. field as 
30 : 9,480 = 0.0045. For 5 fields considered by Kim et 
al., the combined probability for such a fake overdensity 
not to be associated with a QSO is ~ 0.022, and hence 
can be fully neglected by us. Of course they are substan- 
tial over a much larger area considered by Overzier et 
al. 

Alternatively to the low number statistics case, we 
should take into account that Kim et al. found one 2a- 
3(7 overdensity, one — a 2a underdensity field, and three 
additional ones, with no difference to the GOODS fields 
(using highest S/N column in Table 1). If taken at face 
value (that is neglecting likely small number statistics), 
their observational results are consistent with DM halo 
clustering (i.e., high density regions) and an overpopula- 
tion of galaxies, on one hand, and with the scenario in 
which additional effects (i.e., ionization or hydrodynam- 
ics) could play a more dominant role. Accepting results 
of pure DM simulations which show a clear overdensity 
of halos in the QSO field, the problem can lie with our 
understanding of various feedback mechanisms based on 
evolution of starbursts and the central SBHs, i.e., in the 
plausible competition between the positive and negative 
feedbacks which follow from a number of physical pro- 
cesses. This conclusion is weakened by the apparent co- 
incidence between the three QSO fields and the GOODS 
fields (Tabled]). 

The fact that some of the observed QSOs are located 
in average density regions does not necessary imply that 
they are also average DM density regions. For example, 
they could be surrounded by baryon-poor DM halos. A 
variety of negative feedbacks can contribute in removing 
baryons from DM halos or halting the star formation at 
those redshifts, such as the UV background radiation and 
re- ionization (depending on when it occurred) , tidal and 
ram-pressure (ablation) stripping, and dynamical friction 
could play a role. On the other hand, highly and poorly 
collimated outflows from QSOs can provide a positive 
feedback and even trigger the star formation. The inclu- 
sion of baryons in DM simula tions can affect the prop- 
erties of DM structures (e.g., I Romano- Diaz et al.ll2009l 
2010; Du ffy et aI1l2010l ). High resolution numerical sim- 
ulations will be necessary to investigate these and other 
effects. 

Do the Ki m et aTl density estimates reflect the ac- 
tual local densities or "projected" density enhancements 
along the pencil surveys? Observational evidence from 
other QSOs at lower redshifts seems to confirm the Kim 
et al. c onclusions — the QS Os reside in different environ- 
ments (|Stevens et al.|[2010T ). In fact, Kim et al. reports 
a con fidence estimate at the 95% level. iMaselli et al.l 
(2009), using a different method than Kim et al., have ar- 
rived at similar conclusions. Of course, it is not straight- 
forward to relate the mass accumulation and environ- 
ment evolution of high-z QSOs with those at low red- 
shifts. 

In a recent development, Utsumi et al. (2010) have 
used the large field of view of the Suprime-Cam to image 
of the most distant QSO J2329-0301 at z = 6.42 and an 
empty field for comparison. This field of view is about 
two orders of magnitude larger than the HST/ACS field 
of view we have considered here (approximately 0.25 deg 2 



QSOs, the dark side 



9 



versus approximately 11 arcmin 2 ). Seven LBG candi- 
dates have been found in the QSO field compared to only 
one LBG in the comparison field. The authors point out 
that the statistical significance of this apparent overden- 
sity is difficult to estimate. Some evidence exists that the 
LBGs near the QSO avoid the 3 Mpc region centered on 
the QSO — if confirmed this hints to a substantial feed- 
back from the QSO. Extension of this study to additional 
z *~ 6 QSOs would be very useful to confirm widespread 
presence of the galaxy overdensity over large scale. After 
all, the first report on the HST/ACS imaging campaign 
of z ~ 6 QSOs showed a clear overdensity of LBGs for 
the single field analyzed (Stiavelli et al. 2005). 

We now return to the more pessimistic conclusion from 
our comparison between the observational galaxy counts 
and their numerical predictions from our simulations. 
Can this discrepancy be bridged if one assumes bias be- 
tween the DM and baryon distribution? Our assumption 
that there is only one galaxy per DM halo could indeed 
be an oversimplification. However, an increase in the 
number of galaxies per host halo will only aggravate the 
problem, as this will increase the numerical counts. 

An interesting by-product of our analysis of the MAH 
curves for the QSO halos is the observation that the most 
massive halos at high rcdshifts do not necessarily remain 
the most massive ones at lower redshifts, as exhibited by 
Fig. H (see also De Lucia et al 2007; Trenti et al. 2008; 
Overzier et al. 2009). There are corollaries when looking 
for possible low-z counterparts of high-z objects. There 
are a number of issues related to this problem. First, will 
the SBH grow in tandem with its host DM halo? If yes, 
we expect SBH masses of the order of a fewxlO 10 Mq in 
the contemporary universe. Using the comoving volume 
density of high-z QSOs we estimate at least one such 
ultra-massive SBH within a sphere of z <S 0.2. There are 
no observational constraints at present to rule out such 
possibility. 

However, the situation is aggravated because the most 
massive halos today have not been most massive at 
z ~ 6, as our simulations show (see also iTrenti et al.1 
2008 for an analysis of the Millennium simulation merger 
tree history and EPS modeling). In turn, this implies 
that in the local universe there are SBHs more massive 
than the descendants of the z ~ 6 QSO s, if the SBH 
mass is correlated with the D M halo mass (El-Za nt et al.l 
2003||Ferr arese fc Fo rd 2005). That may be contradicted 
by observations as such massive objects will profoundly 
change the galactic dynamics because the radius of influ- 
ence of a fewxlO 10 Mq SBHs will be of the order of the 
visible galaxy. 

6. CONCLUSIONS 

We have performed a set of carefully designed DM 
./V-body simulations aimed at studying the environment 
of QSO-host halos at z ~ 6. A set of 10 Constrained 
Realizations has been constructed, seeding them with a 
10 12 Mq DM halo designed to collapse by z ~ 6 in the 
top-hat scenario. As a control sample, we performed a set 
of 10 simulations that represent the Unconstrained Re- 
alizations of the same CR set. We have constructed the 
halo Mass Functions from each simulation and showed 
how on average, the QSO (i.e., constrained) sample en- 
hanced the DM halo formation in its vicinity due to pure 
gravitational effects, when compared with its respective 



UCR sample. Assuming that there is no bias in baryon 
distribution with respect to the DM in our simulations 
and that the QSO duty cycle is unity, we have calcu- 
lated the expected Lyman Break Galaxies number counts 
and compared th is against the observed QSO fields of 
iKim et all ([200^ . 

Our main result is that the pure DM numerical sim- 
ulations predict a strong o verdensity a round the QSO 
peaks. The observations of IKim et al.l (|2009l) either do 
not support this or provide a mixed result with one pos- 
sibly overdense field, one weakly underdense, and three 
having the density of the GOODS field. The explanation 
for this discrepancy can lie in the small number of QSO 
fields observed coupled with the small number of galaxies 
detected per field. To rule out such possibility it would 
be important to extend the Kim et al. study to a larger 
number of z ~ 6 QSO fields. 

If, however, the over/underdensity result will be con- 
firmed by higher S/N observations, the resolution is prob- 
ably related to the complex physics of feedback mecha- 
nism involving the QSO and possibly the starbursts in 
the surrounding galaxies, or both. With the present anal- 
ysis we are not capable (nor it was our intention) to es- 
tablish the role of feedback in the QSO-host galaxy for- 
mation evolution. However, our results suggest that in 
the case of high density QSO environments, the absence 
of a comparable enhancement of galaxy counts or their 
lack would point to a strong radiative and/or mechanical 
negative feedback from the QSOs, resulting in a strong 
biasing between the DM and baryon distribution. 

A simple analysis of the halo-growth in our CR sample 
together with their MAHs, indicates that these halos will 
possess, by z = 0, a mass in the range of ~ 10 14 M Q . This 
is consistent with our choice of a perturbation mass at 
the initial conditions. 

The underdense field J1148 of E im et al.l that is es- 
pecially at odds with our numerical expectations, opens 
the possibility to future investigations adding baryons to 
simulations to explore the roles of feedback in the evolu- 
tion of QSO environment, star formation, among other 
relevant physical phenomena. 

But perhaps, not all high-z QSOs are formed in high 
density peaks but in less massive (and more common) 
halos of ~ 10 11 Mq. The MAHs of our models show 
that there is natural dispersion (cosmic variance) for the 
halo masses at each z — halos which are most massive at 
z ~ 6 in fact appear least massive at z ~ 3. Hence, QSOs 
could be surrounded by a more "typ ical" enviro nment, 
like in the case of J1048+4637 of the|K im et al.l sample 
(see a similar conclusion by De Lucia & Blaizot (2007) 
and Overzier et al. [2009]). If this is the case, a scenario 
in which negative feedback dominates might not be the 
only solution. Alternatively, z ~ 6 QSOs residing in 
lower mass halos will require evolution of the M. — a 
correlation with z. 
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APPENDIX 

CONSTRAINED REALIZATIONS: NEW FORMALISM 

One of the key ingredients of the canonical model is that the primordial perturbation field constitutes a random 
Gaussian field. Hence the setting of initial conditions (ICs) for cosmological simulations amounts to sampling such 
fields. Gaussian random fields are uniquely determined by their power spectrum, and the assumed statistical homo- 
geneity of the universe, and hence of the perturbation field, implies that Fourier modes of the field are uncorrelated. 
It follows that an unconstrained realization a of Gaussian field can be performed by sampling the independent Fourier 
modes from a normal distribution. Here an unconstrained realization of a Gaussian filed means that the realization 
obeys no other constraint but the assumed power spectrum. 
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A constrained realization (CR) of a Gaussian field is a random realization of such a field constructed to obey a set 
of linear constraints imposed on the field. The imposed constraints violates the statistical homogeneity of th e field , 
hence the random sampling of the Fourier modes cannot be used to construct CRs. The Ho ffman fe; Ribakl (|1991l 
hereafter HR) prov i ded the optim al algorithm for the construction of CRs and it used here for setting the ICs (see 
iZaroubi et allfl995t lHoffmadl2009h . 

The aim is to construct random realizations of a Gaussian field, defined by the WMAP5 power spectrum, constrained 
to have a local density maximum on a mass scale M s centered at the center of the computational box and designed 
to collapse by a redshift z co \\. The li near overdensity co rresponding to the assumed redshift of collapse is calculated 
by the non-linear top-hat model (e.o.. lLokas et al]|2004 and refs. therein). The field under consideration is the linear 
fractional overdensity, <5(r), and the constraints corresponds to the 8 field smoothed with a kernel corresponding to 
M s , A s (r). A local density maximum is defined by 10 constraints: the value of A s , the three components of its 
gradient and the 6 independent components of the Hessian of the field. The three first derivatives are s et to zero, and 
the s econd derivatives are set according the theory of the Gaussian statistics of local density maxima (Bardccn et al. 
1986) . This can be fully integrated into the HR algorithm, but it introduces a considerable computational complexity 
(|van de Wevgaert fe Bertschingerlll996f) . Alternatively the local maximum can be constructed by constraining the 
value of A s at the location of the peak and on six points located along three orthogonal directions, at distances of 
the order of the smoothing length corresponding to M s . Given the three angles defining the three orthogonal axes 
th ese constitute 10 con straints. The off-peak values can be calculated by using the full Gaussian statistics machinery 
of iBardeen etafl (|1986l ). 

Before the proceeding to the construction of CRs the smoothing kernel needs to be specified. Analytic considerations 
favor the use of a top-hat filter. However, the sharp real space cut-off leads to the Gibbs Phenomenon, the so-called 
in k-space, namely the 'ringing' of corresponding fields in k-space. To avoid that a Gaussian filter is adopted here. 
Namely, 

A s « = * / cxpf-^-# W)dV. (Al) 



(27Ti?2)3/2 J 2i? 2 

The smoothing length is related to M s by 

r -= om {W- <A2) 

where p is the mean cosmological density ()Bardeen et al.|[l~986l) . 

Given a random realization 8(r) and a set of constraints {A s (r Q ,) ct=lj 7 }, where r a is the location at which the a 
constraint is imposed, a CR is obtained by: 

5 CR (r)=5(r)+aA 3 (r-r Q )(eA s A 3 (rar^) + a^ /3 )^ 1 (A s (^)-A s (r^)), (A3) 



where P(k) is the power spectrum. The various terms that appear in Eq |A3l are defined as follows. The cross-correlation 
function of the underlying and the smoothed 8 fields is given by: 

&,A 8 = \h- [ exp(-^-)p(k)k 2 dk. (A4) 



The auto-correlation function of the smoothed field is given by: 

£a s ,a 3 = \ - / cxp ( - k 2 R 2 s )P(k)k 2 dk (A5) 



The pseudo constraints A s (r^) are obtained by convolving the random S(r) with the Gaussian filter, evaluated at the 
position of the constraints. The <J 2 a 8^ a term in the inverse of the auto-correlation matrix adds the possibility of a a a 
'fuzziness', or uncertainty, in the value of the imposed constraint A s (r Q ). (6^ p is the Kronecker 8.) 

The mean density field around a given density maximum is determined by the assumed power spectrum, the height 
of the maximum and the Hessian matrix (Bardcc n et al.lll986D . It can be shown that for high peaks the mean density 
profile around the peak is closely approximated by the mean field around a random field point. This property is used 
here to calculate the values of the off-peak constrai nts. The mean field around a filed point, constrained to obey just 
A s (r = 0), is given by the Wiener Filter estimator (|Zaroubi et al.l lT995; Hoffman 2009), which is easily obtained from 
Eq. lA3l bv setting 8(r) to zero. It follows that the construction of a CR proceeds in two steps. First, a single constraint 
is imposed and the mean field given the constraint is constructed. The values of the off-peak constraints are extracted 
from the mean field. Then, having the full set of constraints at hand the ensemble of CRs is constructed. 

The present impl e mentation of the CR algorithm improves on the earlier utilization of the method employed by 
iRomano-Diaz et al.l (|2006l . 120071 120081) . There DM halos were imposed by using one single constraint for a given halo, 
i.e. halos have been treated as random field points as far as the constraints are considered. The present implementation 
results in better constrained halos, which faithfully obey the imposed (linear) constraints. In particular the virial mass 
of the imposed halo is very close to the imposed M s and it is situated very close to the position imposed by the 
constraints. 
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# X [Mpc fe" 1 ] Y [Mpc h- 1 ] Z [Mpc A s a 



1 


0.0 


0.0 


0.0 


9.05 0.0 


2 


-0.94 


0.0 


0.0 


8.39 0.1 


3 


0.94 


0.0 


0.0 


8.39 0.1 


1 


0.0 


-0.94 


0.0 


8.39 0.1 


5 


0.0 


0.94 


0.0 


8.39 0.1 


6 


0.0 


0.0 


-0.94 


8.39 0.1 


7 


0.0 


0.0 


0.94 


8.39 0.1 



TABLE Al 

The positions and values of the constraints used to construct the set the initial conditions of the 

simulations. 



The set of constraints used for a mass scale M s — 10 12 M Q /i _1 , where h is Hubble's constant in units of 
100 km s -1 Mpc -1 , collapse redshift of z co \\ = 6.0 and the assumed WMAP5 cosmology is presented in Table lATl For 
simplicity the realizations are evaluated at the present epoch, z — 0, and are then scaled down their initial values by 
the linear growth factor of the WMAP5 cosmology. 



